February 7, 2008 23:33 



WSPC - Proceedings Trim Size: 9in x 6in cortona06 



1 



Quantum Monte Carlo calculations of symmetric nuclear matter 

S. GANDOLFI and F. PEDERIVA 

Dipartimento di Fisica and INFN, University of Trento, via Sommarive 14, 1-38050 

Povo, Trento Italy 

S. FANTONI 

S.I.S.S.A., International School of Advanced Studies and INFN via Beirut 2/4, 34014 

Trieste, Italy 

K. E. Schmidt 

Department of Physics, Arizona State University, Terape, AZ, USA 



We present an accurate numerical study of the equation of state of nuclear mat- 
ter based on realistic nucleon— nucleon interactions by means of Auxiliary Field 
Diffusion Monte Carlo (AFDMC) calculations. The AFDMC method samples 
the spin and isospin degrees of freedom allowing for quantum simulations of 
large nucleonic systems and represents an important step forward towards a 
quantitative understanding of problems in nuclear structure and astrophysics. 

1. INTRODUCTION 

The knowledge of the properties of symmetric and asymmetric nuclear mat- 
ter, is one of the main ingredients of the equations predicting the structure, 
the dynamics and the evolution of stars, in particular during their last 
stages, when they become ultra-dense neutron stars. The large variety of 
available results, partly due to the use of many different approximations, 
and partly to the limited knowledge of the nuclcon-nuclcon interaction, 
lead to different interpretations and hypothesis. Depending on the model 
used, the density of nuclear matter in the inner shells can reach up to 9 
times the core density of stable nuclei, po = 0-16 fm -3 , 1 or be limited to 
4-5 times po- 

One important step towards the understanding of these astrophysical 
problems is the study of the symmetric nuclear matter, related to the var- 
ious model NN interactions available. While considerable advances have 
been made, 2:3 it is still impossible to understand how different results de- 
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pend on the model interaction used rather than on the particular approxi- 
mation employed for computing energies. 

Therefore, it is not possible to draw conclusions when comparing with 
astronomical observations. However, the recent success in predicting the 
properties of light nuclei gives us some confidence that the non-relativistic 
description of nuclear matter based on effective potentials fitted to repro- 
duce NN data and the binding energy of light nuclei can be reliable enough. 
The main feature of such nucleon-nucleon interactions, besides the short 
range repulsion, is the explicit dependence on the relative quantum state 
of the nucleons which can be described using spin and isospin, angular 
momentum, spin-orbit, and tensor operators. 4 

Properties of light nuclei (A<6) can be efficiently computed with high 
accuracy using modern few-body techniques 5-7 or with the ab initio no-core 
nuclear shell model (with A<12 8 ). Quantum Monte Carlo techniques based 
on recasting the Schrocdinger equation into a diffusion equation (Diffusion 
Monte Carlo or Green's Function Monte Carlo), allowed for performing 
calculations up to A<12. 9 ' 10 However, the computational resources needed 
for such simulations are very large, because of the summation over all the 
possible states necessary to evaluate the terms of the Hamiltonian with 
a quadratic dependence on spin and isospin. The number of such terms, 
and the CPU time needed to calculate them, grows exponentially with the 
number of nucleons; 12 C 9 or 14 neutrons 11 is the limit for the currently 
available computational resources. 

Since the spatial degrees of freedom are already sampled, one would like 
to replace the sum over the spin-isospin states with an efficient sampling 
method. The simulation of symmetric nuclear matter requires a minimum 
of 28 nucleons in a box replicated in space (7 nucleons for each spin-isospin 
state) to obtain a wave function with closed shells of momenta, and this is 
out of the reach of the standard Quantum Monte Carlo methods. 

Recently 12 we shown that the spin-isospin is efficiently sampled by using 
the Auxiliary Field Diffusion Monte Carlo method, 13 which is based on the 
use of auxiliary variables to linearize the quadratic spin-isospin operators of 
the nuclear matter Hamiltonian, making them treatable in a diffusion Monte 
Carlo scheme. Up to now it has been applied to simulate pure neutron 
matter (up to 114 neutrons), 14,15 and neutron drops 16 ' 17 interacting with 
realistic two- plus three-body interactions. 

Here we extend calculations to include isospin degrees of freedom and to 
deal with the strong tensor-isospin force, responsible of the nuclear binding. 
The method can readily handle an asymmetry in the number of neutron 
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and protons or the deformation of heavy nuclei. 
2. METHOD 

Auxiliary Field Diffusion Monte Carlo (AFDMC) 13 is an extension of 
the standard Diffusion Monte Carlo method in which the ground state 
of an Hamiltonian H is obtained by solving the imaginary time depen- 
dent Schroedinger equation. The solution is obtained by evolving a pop- 
ulation of configurations of the system ("walkers") X = {R,a,r}, where 
R = {fi, ...,f N ,},a = {f?i, . . . ,<tjv}, and r = {n, . . .,t n }, with F(X,t) = 
^ T (X)^(X,t), according to 

F(X,t) = J dX'^^G (X,X',t)F(X',0) (1) 

The function is a "trial" wave function, usually determined by means 
of variational calculations, and Go is an approximation to the Green's func- 
tion of the imaginary time Schroedinger equation: 

G (X,X',t) = ^ Dt )-^/2 e -(R-R'f/ADt e -t{V(X)-E )^ (2) 

where D = h 2 /2m, E is an estimate of the ground state energy of the 
system, and V(X) is the nucleon- nucleon interaction. For a long enough 
imaginary time the distribution of the walkers converges to the product 
^>T(X)^fo(X) where is the wave function of the ground state of H. This 
fact allows the computation of matrix elements (<I't|0| , I'o) of observables 
O of interest in a Monte Carlo way. When O = H the value obtained is the 
exact ground state energy of the system. 

The presence of an interaction V(X) including operators like (3(7^ ''f'ij^j ' 
fij — <Ji ■ <Jj ) and Tj • fj is the origin of the computational cost in the standard 
approaches. AFDMC uses the Hubbard-Stratonovich method to transform 
the operators S which are quadratic in the spin and isospin into linear 
operators: 

e -(i/2)txs- = 1 / dye -v 2 /\vV^is (3) 

Then S are operators which are linear combinations of the spin and 
isospin operators for each nucleon, and A depend on the interaction. The 
transformed Green's Function is applied to the spin-isospin part of the wave 
function, and its effect consists of a rotation of the spin and isospin degrees 
of freedom (written as four-component spinors in the proton-neutron up- 
down basis) by a quantity that depends on the auxiliary variable y along 
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with multiplication of the state by an overall factor. The sum over spin and 
isospin is replaced by sampling a set of rotations of the variables. This pro- 
cedure reduces the dependence of the computational time on the number 
of nucleons necessary for performing a simulation step from exponential to 
cubic. It is therefore possible to perform on a regular workstation or on a 
modest PC cluster calculations that would require Tflop supercomputers 
with the standard methods. This method, like other diffusion Monte Carlo 
methods, suffers from the so-called "sign problem" when it is applied to 
fermions, and when complex wave functions need to be used. In our calcula- 
tions we apply the fixed-phase approximation to overcome this problem. 18 
While this method has already been successfully applied to pure neutron 
matter, 15 it has not been previously used for mixed proton and neutron 
systems. It should be noted that it does not guarantee an upperbound to 
the mixed energy used here. As a test for the correctness and the efficiency 
of our approach we reproduced within 0.3 MeV total energy the existing 
results for the binding energy of 4 He with potentials of the AV6 class. 19 
We have also been able to compute binding energies for 16 O and for 40 Ca 20 
with this method. 

3. INTERACTION 

A crucial point in dealing with nuclear matter is the choice of the inter- 
action among nucleons. As already mentioned, several modern two-body 
potentials are available nowadays, all fitting the NN data with x 2 ~ 1- We 
use the potentials of the Argonne class with n operators (AVn). 21 While 
the full version contains n = 18 operators, most of the physics is reasonably 
well described by the first 6 operators made up of 4 central spin-isospin 
dependent components and two tensor ones, which include the long range 
one-pion-exchange force. The most important missing terms are the spin- 
orbit components. In nuclei and neutron drops the spin-orbit contribution 
to the energy amounts to a few tenths of MeV/nucleon. A correction of the 
order 1 MeV/nucleon can be attributed at low densities to the remaining 
terms included in AV18. Specifically, we have used the interaction Argonne 
i>g 19 truncated by dropping the spin-orbit terms, and including only the 
first six operators, which we denote as AV6'. 

It is well known that two-body NN interaction underbind light nuclei, 
and one needs to add a specific effective three-body potential to reproduce 
their low energy properties. Semi-phenomenological three-nucleon interac- 
tions following the lowest order three-nucleon diagrams with one and two 
intermediate Delta resonance states provide a very satisfactory description 
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of the ground state energy and the low level spectra of light nuclei up to 
12 C. 9 We have disregarded such three-body forces in our simulations. In 
nuclear matter they are essential to reproduce the experimental saturation 
density, and, in general, they contribute about 10% of the total binding 
energy. A full comparison with the available experimental data goes be- 
yond the scopes of the present paper. Here we are interested in showing the 
efficiency of the AFDMC methods in dealing with nuclear matter models 
which include realistic tensor interactions like in the AV6' potential. 

4. RESULTS 

The results of the calculations with A=28 include box corrections that have 
computed by adding to the two body sums contribution of nucleons in the 
first shell of periodic cells. Such procedure is effective. In order to assess 
the magnitude of finite size effects we performed calculations with 76 and 
108 nucleons at densities p = 0.08 fm~ 3 and p = 0.48 fm~ 3 . The results 12 
coincide with the ones obtained with 28 nucleons within 3 percent. 

We computed the EOS of symmetric nuclear matter in the range of 
densities 0.5 < (p/po) < 3, and compared it with previous available results 
obtained with the same potential using Fermi Hypernetted Chain in the 
Single Operator Chain approximation (FHNC/SOC) and the Brueckner- 
Hartrcc-Fock (BHF) in the two-hole line approximation. 22 AFDMC cal- 
culations were performed with 28 nucleons, filling the shell of plane waves 
with momentum of modulus 1 and providing a wave function yielding an 
isotropic density. 

The results are summarized in Fig. 1 and reported in Ref. 12. The com- 
parison of the various EOS suggests the following comments: FHNC/SOC 
leads to an ovcrbinding at high density. A similar indication was found by 
Moroni et al. 23 after a DMC calculation of the EOS of normal liquid 3 He 
at zero temperature, with a guiding function including triplet and backflow 
correlations. The comparison with the equivalent FHNC/SOC calculations 
of Refs. 24,25 have shown similar discrepancies. There are two main intrin- 
sic approximations in variational FHNC/SOC calculations, which violate 
the variational principle. The first one consists in neglecting a whole class 
of cluster diagrams, the so called elementary diagrams, which cannot be 
summed up by means of FHNC integral equations. We have calculated the 
lowest order diagram of this class, namely the one having only one correla- 
tion bond and four exchange bonds. The results obtained show a substantial 
effect from this diagram and bring the FHNC/SOC estimates very close to 
AFDMC results, as shown in Fig. 1. The second approximation is related 
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Fig. 1. (color online). Equation of state of symmetric nuclear matter calculated with 
different methods. Red circles represent AFDMC results with statistical error bars and 
the green line is the fitted functional form described in the text. Dashed lines correspond 
to calculations performed with other methods 22 (blue line with squares: FHNC/SOC; 
magenta with diamonds: BHF). Blue triangles represent the FHNC/SOC energies cor- 
rected by including the low order of elementary diagrams as described in the text. Blue 
arrows show the corresponding energy shift, which increases at higher densities. 



to the non-commutativity of the correlation operators entering the varia- 
tional wave function. The only class of cluster diagrams contributing to such 
non-commuting terms, which can be realistically calculated, is that charac- 
terized by single operator chains. It is believed that such an approximation 
is reliable in nuclear matter, but there is no clear proof of this. 

BHF calculations of ref. 22 predict an EOS with a shallower binding than 
the AFDMC one. It has been shown for symmetric nuclear matter, using 
the AV18 and AV14 potentials, that contributions from three hole-line di- 
agrams add a repulsive contribution up to ~ 3MeV at densities below po, 2e 
and decrease the energy at high densities. 27 Such corrections, if computed 
with the AV6' potential, would probably preserve the same general behav- 
ior, and bring the BHF EOS closer to the AFDMC one. Therefore, our 
calculations show that the two hole-line approximation used in Ref. 22 is 
too poor, particularly at high density. 

The AFDMC equation of state was fitted with the following functional 
form: 

^ = ^+a(x-x) 2 +{3(x-x)\ (4) 
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where x = p/po and the various coefficients are given by Eq/A = -14.04(4) 
MeV, a = 3.09(6) MeV, = -0.44(8) MeV, and x = 1.83(1). The resulting 
compressibility K = 9x 2 (9 2 (E/A) /dx 2 )_ at saturation density x is ~ 190 
MeV. The fit of the EOS allows for computing the pressure vs. density for 
symmetric nuclear matter. 

5. CONCLUSIONS 

The availability of an efficient and relatively fast projection algorithm for 
the computation of energies and other observables of dense hadronic matter 
enables the possibility of a more quantitative understanding of the prop- 
erties of neutron stars and supernovae, as well as that of medium-heavy 
nuclei. Computations on such systems are at present out of reach of the 
standard GFMC methods and available supercomputers. Therefore, the 
extension of AFDMC algorithm to deal with nuclear matter is a significant 
step forward. Some technical improvements on the calculations presented 
here, such as the addition to the AV6' of spin-orbit terms and three-body 
interactions are already underway. The treatment of asymmetric nuclear 
matter, particularly important for the determination of the properties of 
neutron stars, is also straightforward, and will be the subject of future 
exploration. 

We acknowledge helpful conversations with M.H. Kalos, P. Faccioli, W. 
Leidemann, E. Lipparini, and G. Orlandini. This work was in part sup- 
ported by NSF grant PHY-0456609. Calculations were performed on the 
HPC facility "BEN" at ECT* in Trento under a grant for Supercomputing 
Projects. 
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